clear all
numIter = 10;
word_len = 7;
period_labels = [10 20 30];
numTargets = 3;


best_baseSeq = [];
best_probability_repetition = [];
best_transition_sequence = [];
best_strSeq = [];
best_fullSeq = [];
best_transSeq = [];
best_detectPower = [];

ITIlist = 1500:100:10000;

for j=1:length(ITIlist)
    
    ITI = ITIlist(j);
    maxDetectPower = 0;

    for i=1:numIter
        fprintf('Iteration %d, ITI=%d:\n',i, ITI);
        % Generate sequence using a sinusoidal guide function with period=period_labels
        [baseSeq probability_repetition transition_sequence] = seqGen(word_len,period_labels);
        
        % Add zeros and targets
        switch numTargets
            case 3
                [strSeq, fullSeq] = addZerosAnd3Targets(baseSeq);
            case 4
                [strSeq, fullSeq] = addZerosAnd4Targets(baseSeq);
            case 5
                [strSeq, fullSeq] = addZerosAnd5Targets(baseSeq);
        end
        
        fullSeq = fullSeq(:);
        
        % Generate transition sequence
        transSeq = zeros(size(fullSeq));
        transSeq(1) = -1; % The first transition is invalid (blank to blank)
        for j=2:length(transSeq)
            if min(fullSeq(j),fullSeq(j-1)) < 1
                transSeq(j) = -1;
            else
                transSeq(j) = abs(fullSeq(j)-fullSeq(j-1));
            end
        end
        
        % Calculate detection power
        detectPower = evalSeq(transSeq, ITI, 'hrf.txt', 10);
        
        if detectPower > maxDetectPower
            %{
        display(baseSeq');
        display(probability_repetition);
        display(transition_sequence);
        display(strSeq);
        display(fullSeq);
        display(transSeq);
            %}
            
            best_baseSeq = [best_baseSeq baseSeq];
            best_probability_repetition = [best_probability_repetition probability_repetition'];
            best_transition_sequence = [best_transition_sequence transition_sequence'];
            best_strSeq = [best_strSeq strSeq(:)];
            best_fullSeq = [best_fullSeq fullSeq];
            best_transSeq = [best_transSeq transSeq];
            best_detectPower = [best_detectPower detectPower];
            
            display(detectPower)
            if detectPower > maxDetectPower
                maxDetectPower = detectPower;
            end
        end
        
    end

end


save all